#> Info: preparing to read 1 data file(s)...
#> Info: reading file 1/1 'Garrett_Adducts.cf.rda' with '.rda' reader...
#> Info: loaded data for 71 data files from R Data Archive - checking loaded files for content consistency...
#> Info: encountered 6 problems in total.
#> # A tibble: 6 x 4
#> file_id type func details
#> <chr> <chr> <chr> <chr>
#> 1 587__Linearity CO2__.dxf error extract_dxf_r… cannot identify measured …
#> 2 587__Linearity CO2__.dxf error extract_isoda… `peak`, `start`, `rt`, `e…
#> 3 586__Linearity CO2__.dxf error extract_dxf_r… cannot identify measured …
#> 4 586__Linearity CO2__.dxf error extract_isoda… `peak`, `start`, `rt`, `e…
#> 5 585__Linearity CO2__.dxf error extract_dxf_r… cannot identify measured …
#> 6 585__Linearity CO2__.dxf error extract_isoda… `peak`, `start`, `rt`, `e…
#> Info: aggregating vendor data table without units from 71 data file(s), including file info 'c(ox = `Seed Oxidation`, everything())'
# add metadata
data_w_metadata <- vdt %>%
iso_add_metadata(metadata_samples, match_by = c(`Identifier 1`))
#> Info: metadata added to 1827 data rows, 40 left without metadata:
#> - 25 metadata entries were mapped to 1827 data entries via column 'Identifier 1'
# missing metadata
data_w_metadata %>%
iso_get_missing_metadata(select = c(Analysis, `Identifier 1`, map_id))
#> Info: fetching data entries that are missing metadata
# map peaks
data_w_peaks_all <- data_w_metadata %>%
# only continue with records that have metadata and are flagged for processing
iso_remove_missing_metadata() %>%
filter(process == "yes") %>%
# map peaks
iso_map_peaks(metadata_peak_maps, file_id = file_id, rt = Rt, rt_start = Start, rt_end = End)
#> Info: removing 40 of 1867 data entries because of missing metadata
#>
# problematic peaks
prob <- data_w_peaks_all %>%
iso_get_problematic_peaks(select = c(file_id, Analysis, compound, Start, Rt, End))
#> Info: fetching 596 of 2044 data table entries with problematic peak identifications (unidentified, missing or ambiguous)
# focus on non problematic peaks only?
data_w_peaks <- data_w_peaks_all %>% iso_remove_problematic_peaks()
#> Info: removing 596 of 2044 data table entries because of problematic peak identifications (unidentified, missing or ambiguous)
data_w_peaks
data_w_peaks %>%
iso_plot_ref_peaks(
is_ref_condition = is_ref_peak == "yes",
ratio = c(`R 13C/12C`, `R 18O/16O`),
group_id = `Analysis`,
within_group = TRUE
)
Focus on the analytes and calculate a few summary parameters we want to use later.
Add known isotope values for standards.
#> Info: matching standards by 'type' and 'compound' - added 15 standard entries to 458 out of 762 rows
Determine calibrations fits for individual standard analyses.
#> Info: preparing data for calibration by grouping based on 'Analysis'
#> Info: generating 'd13C' calibration based on 1 model ('lm(`d 13C/12C` ~ true_d13C)') for 31 data group(s) with standards filter 'is_standard'. Storing residuals in new column 'd13C_resid'.
#> Info: there were no problematic calibrations
#> Warning in `[<-.data.frame`(`*tmp*`, is_list, value = structure(list(`3` =
#> "<>", : replacement element 1 has 1 row to replace 0 rows
#> Warning in `[<-.data.frame`(`*tmp*`, is_list, value = structure(list(`3` =
#> "<>", : replacement element 2 has 1 row to replace 0 rows
Look at calibration parameters (coefficients and summary) as well ad data ranges in overview table.
Visualize the calibration parameters
# unnest some of the data for plotting
data_for_plots <-
stds_w_calibs %>%
iso_unnest_data(
c(Oxidation = `Seed Oxidation`,
datetime = file_datetime,
`Mean Amplitude` = ampl_sample_mean.mV)
)
data_for_plots
# parameters vs time
data_for_plots %>%
iso_plot_calibration_parameters(
calibration = "d13C",
x = datetime,
color = `Injection Volume`,
shape = Oxidation
)
# parameter vs analysis
data_for_plots %>%
iso_plot_calibration_parameters(
calibration = "d13C",
x = paste(Analysis),
color = `Injection Volume`,
shape = Oxidation
)
# parameters vs amplitude
data_for_plots %>%
iso_plot_calibration_parameters(
calibration = "d13C",
x = `Mean Amplitude`,
color = `Injection Volume`
)
# turn the last plot into an interactive one
ggplotly(ggplot2::last_plot() + theme(legend.position = "none"))
Look at residuals to see goodness of fit for the single analysis calibrations.
stds_w_calibs %>%
# pull out all data
iso_unnest_data(select = everything()) %>%
filter(is_standard) %>%
# calculate deviation from means
group_by(Analysis) %>%
mutate(
`Var: residual d13C [permil]` = d13C_resid,
`Var: area diff from mean [%]` = (`Intensity 44`/mean(`Intensity 44`) - 1) * 100,
`Var: amplitude diff from mean [%]` = (`Ampl 44`/mean(`Ampl 44`) - 1) * 100
) %>%
ungroup() %>%
# visualize
iso_plot_data(
x = compound,
y = starts_with("Var"),
group = Analysis,
color = `Injection Volume`
) +
# plot modifications
labs(x = "")
ggplotly()
Run global calibrations across all standards.
#> Info: preparing data for calibration by nesting the entire dataset
#> Info: generating 'd13C' calibration based on 4 models ('delta_only = lm(`d 13C/12C` ~ true_d13C)', 'delta_and_ampl = lm(`d 13C/12C` ~ true_d13C + `Ampl 44`)', 'delta_cross_ampl = lm(`d 13C/12C` ~ true_d13C * `Ampl 44`)', 'delta_and_time = lm(`d 13C/12C` ~ true_d13C + file_datetime)') for 1 data group(s) with standards filter 'type == "standard" & is_standard'. Storing residuals in new column 'd13C_resid'.
Visualize the calibration parameters for the different models.
# coefficient summary table
data_w_calibs %>% iso_unnest_calibration_parameters("d13C")
# plot to look at coefficients and summary
data_w_calibs %>%
iso_plot_calibration_parameters(
calibration = "d13C",
x = d13C_calib,
color = d13C_calib,
# include the significance level to gauge which model parameters matter
shape = signif
)
Garrett thinks delta and aplitude and delta cross amplitude look better
data_w_calibs %>%
iso_unnest_data(select = everything()) %>%
rename(`Injection Volume` = `Identifier 2`) %>%
filter(type == "standard" & is_standard) %>%
# plot each standard analysis for each of the models to see the differences
iso_plot_data(
x = compound,
y = d13C_resid,
group = paste(Analysis, d13C_calib),
color = d13C_calib
)
Apply the delta_and_ampl calibration.
#> Info: applying 'd13C' calibration to infer 'true_d13C' for 1 data group(s); storing resulting value in new column 'true_d13C_pred' and estimated error in new column 'true_d13C_pred_se'. This may take a moment... finished.
The comparison of the predicted value in the standards along with the standard error in the prediction (the inherent error of the global calibration itself, which should be treated as the analytical error) can give a better sense of how precise the data is.
Keep in mind that the calibration provides not just the predicted value of a peak but also the standard error based on the calibration and an indicator whether a predicted data point was in the calibration range or not.
data_calibrated %>%
iso_plot_data(
x = Analysis,
y = true_d13C_pred,
y_error = true_d13C_pred_se,
color = true_d13C_pred_in_range,
shape = type,
points = TRUE,
lines = FALSE
) +
# additional features beyond the default plot
facet_wrap(~compound, scales = "free_y") +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(color = "in calib. range")
final_data <- data_calibrated %>%
#filter out bad chromatograms - all are the lower amount injection, which were followed by a larger volume - keeping the larger volume injection data
filter(Analysis != "597") %>% # 122 depth, OG005
filter(Analysis != "593") %>% # 121.105, OG047
filter(Analysis != "612") %>% # 121.89, OG040
filter(Analysis != "613") %>% # 120.205, OG042-2
filter(Analysis != "609") %>% # 119.55, OG158
filter(Analysis != "595") %>% # 122.9, OG043
filter(type == "sample") %>%
filter(true_d13C_pred_in_range == TRUE) %>%
select(
# sample info
depth = `Depth`, `compound`, Identifier.1 = `Identifier 1`, `ampl_sample_mean.mV`, type, Analysis,
# peak info
Ampl.44 = `Ampl 44`, true_d13C_pred , true_d13C_pred_se, true_d13C_pred_in_range, d13C_resid
) %>%
mutate(
prep = "adduct"
)
vis <- final_data %>%
filter(compound %in% c("C29", "C27", "C31", "C33"))%>%
iso_plot_data(
x = depth,
y = true_d13C_pred,
y_error = true_d13C_pred_se,
size = Ampl.44
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
scale_x_reverse() +
coord_flip()
ggplotly(vis)
lowres <- read.csv("lowres_d13c.csv") %>% select( -X, -`measurement_info`, -`file_datetime`, -file_id) %>% mutate(prep = "non-adduct")
sapply(final_data, class)
#> depth compound Identifier.1
#> "numeric" "character" "character"
#> ampl_sample_mean.mV type Analysis
#> "numeric" "character" "character"
#> Ampl.44 true_d13C_pred true_d13C_pred_se
#> "numeric" "numeric" "numeric"
#> true_d13C_pred_in_range d13C_resid prep
#> "logical" "numeric" "character"
sapply(lowres, class)
#> compound Identifier.1 depth
#> "factor" "factor" "numeric"
#> ampl_sample_mean.mV type Ampl.44
#> "numeric" "factor" "numeric"
#> true_d13C true_d13C_pred true_d13C_pred_se
#> "numeric" "numeric" "numeric"
#> true_d13C_pred_in_range d13C_resid prep
#> "logical" "numeric" "character"
final_data <- bind_rows(lowres, final_data)
#> Warning in bind_rows_(x, .id): binding factor and character vector,
#> coercing into character vector
#> Warning in bind_rows_(x, .id): binding character and factor vector,
#> coercing into character vector
#> Warning in bind_rows_(x, .id): binding factor and character vector,
#> coercing into character vector
#> Warning in bind_rows_(x, .id): binding character and factor vector,
#> coercing into character vector
#> Warning in bind_rows_(x, .id): binding factor and character vector,
#> coercing into character vector
#> Warning in bind_rows_(x, .id): binding character and factor vector,
#> coercing into character vector
final_data <- final_data %>%
filter(type == "sample") %>%
select(-true_d13C, -d13C_resid) %>%
filter(compound %in% c("phytane", "C17", "C19", "C21", "C27", "C29", "C31", "C33", "C35")) %>%
# net fractionation factor for inorganic -> lipid d13C (eps_net) calculated as = (inorganic -> biomass fractionation factor) + (biomass -> lipid fractionation factor)
# minimum and max for marine algal biomass (first number in equation) from spread in data of algal biomass in Ohkouchi et al., 2015; OPTIONS FOR SHORTCHAIN: -4 per mil for phytane (Schouten et al., 1998; Sinninghe Damste et al., 2008), -8.4 per mil for cyano nC17 , -6.4 per mil for phytane (above takes into account frctionation in sulfur-bound and release per Sinninghe Damste?)(Sakata et al., 1997)
# lipid-leaf addition (second number in equation) by chain length from gymnosperm data in Diefendorf et al., 2011 [NOTE: average for gymnosperm is -2.5 per mil, that may be more appropriate]; minimum and maximum offset of leaf to atmosphere from C3 plant histogram in Ehleringer and Cerling, 2001 (in Taiz et al., 6th ed. Plant Physiology and Development)
mutate(eps_net_min = case_when(.$compound %in% c("C17", "C19", "C21") ~ (12+4) ,
.$compound == "phytane" ~ (12+4.4),
.$compound == "C27" ~ (22+2),
.$compound == "C29" ~ (22+3),
.$compound == "C31" ~ (22+4),
.$compound == "C33" ~ (22+3.5),
.$compound == "C35" ~ (22+3)),
eps_net_max = case_when(.$compound %in% c("C17", "C19", "C21") ~ (25+8) ,
.$compound == "phytane" ~ (25+6.4),
.$compound == "C27" ~ (32+2),
.$compound == "C29" ~ (32+3),
.$compound == "C31" ~ (32+4),
.$compound == "C33" ~ (32+3.5),
.$compound == "C35" ~ (32+3)),
d13C_pool_min = true_d13C_pred + eps_net_min,
d13C_pool_max = true_d13C_pred + eps_net_max,
d13C_pool_avg = true_d13C_pred + ((eps_net_max + eps_net_min) / 2) )
final_data <- final_data %>%
mutate(
d13C_p = case_when(.$compound %in% c("C17", "C19", "C21") ~ true_d13C_pred + 8.4,
.$compound == "phytane" ~ true_d13C_pred + 6.4 ) , # from Julio's unpublished text, original from Sakata et al., 1997 (6.4 in Sakata...)
eps_p_lipid = 1000 * ((d13C_pool_avg + 1000)/ (d13C_p + 1000) -1), # from Freeman and Hayes 1992, but using avg estimated d13C of DIC from lipid rather than CaCO3 measurement. This allows moving pool value - important when pool itself is experiencing characterizable d13C changes thru time.
eps_f = case_when(.$compound %in% c("C17", "C19", "C21") ~ 24,
.$compound == "phytane" ~ 25 ) , # from Julio's unpublished text, 22 for alkanes (cyano-plank intermediate) and 25 for phytane; Pagani 2002 says 25 for phytoplankton in general; "Goericke et al. [ 1994] calculate the likely range of eps_f for marine phytoplankton to be 25-28" (from Brigidare et al., 1997). After Modelling (below), 24 for short chains gets us closest to phytane value (assuming 25 for phytane) - has implications for makeup of community, i.e., that cyanos are lower abundance than phytoplankton
b = 170, # per Julio's unpublished text, can sequence from [118 to 262], or use constant 170. I think at least 170, potentially higher, given modelled need for high [PO4] advection to sustain OAE2
CO2_aq_conc = b / (eps_f - eps_p_lipid) , # from Brigidare et al. 1997
Ko = (2.418 * 10^-2), # solubility constant, from Weiss 1974 (tables II and III); dep. on T and Salinity, range from ~2 to ~7 over all regimes; @ ~35 C, range is from 2.1 to 2.3 across all salinities. Choosing for 34ppt (Hay et al., 1996) and 32C (mid range T from O'Brient et al., 2017)
pCO2 = CO2_aq_conc / Ko # from Julio's unpublished text , and Sarmiento and Gruber, Ocean Biogeochemical Dynamics, section 8.3 (p.327)
)
final_data %>% filter(compound == "C21") %>% select( compound, eps_p_lipid, d13C_pool_avg, true_d13C_pred, CO2_aq_conc, pCO2) %>% kable()
| compound | eps_p_lipid | d13C_pool_avg | true_d13C_pred | CO2_aq_conc | pCO2 |
|---|---|---|---|---|---|
| C21 | 16.45642 | -5.558408 | -30.05841 | 22.53572 | 931.9983 |
| C21 | 16.44228 | -4.716870 | -29.21687 | 22.49355 | 930.2542 |
| C21 | 16.44119 | -4.652097 | -29.15210 | 22.49031 | 930.1204 |
| C21 | 16.41665 | -3.188455 | -27.68846 | 22.41754 | 927.1108 |
| C21 | 16.42721 | -3.818942 | -28.31894 | 22.44880 | 928.4037 |
| C21 | 16.42893 | -3.921605 | -28.42161 | 22.45390 | 928.6147 |
| C21 | 16.39894 | -2.129421 | -26.62942 | 22.36531 | 924.9508 |
| C21 | 16.47220 | -6.495370 | -30.99537 | 22.58295 | 933.9514 |
| C21 | 16.45913 | -5.719790 | -30.21979 | 22.54383 | 932.3338 |
| C21 | 16.45878 | -5.698920 | -30.19892 | 22.54278 | 932.2904 |
| C21 | 16.43330 | -4.182080 | -28.68208 | 22.46687 | 929.1508 |
| C21 | 16.45231 | -5.313985 | -29.81398 | 22.52345 | 931.4907 |
| C21 | 16.39212 | -1.721012 | -26.22101 | 22.34527 | 924.1218 |
| C21 | 16.42510 | -3.692964 | -28.19296 | 22.44254 | 928.1449 |
| C21 | 16.39335 | -1.794682 | -26.29468 | 22.34888 | 924.2712 |
| C21 | 16.40776 | -2.656710 | -27.15671 | 22.39127 | 926.0244 |
| C21 | 16.43223 | -4.118406 | -28.61841 | 22.46370 | 929.0197 |
| C21 | 16.48668 | -7.354039 | -31.85404 | 22.62648 | 935.7518 |
| C21 | 16.44045 | -4.608185 | -29.10818 | 22.48812 | 930.0296 |
| C21 | 16.45807 | -5.656454 | -30.15645 | 22.54065 | 932.2021 |
| C21 | 16.44026 | -4.596741 | -29.09674 | 22.48755 | 930.0060 |
| C21 | 16.44202 | -4.701641 | -29.20164 | 22.49279 | 930.2227 |
| C21 | 16.44346 | -4.787331 | -29.28733 | 22.49707 | 930.3999 |
| C21 | 16.43946 | -4.548936 | -29.04894 | 22.48516 | 929.9073 |
| C21 | 16.38667 | -1.394377 | -25.89438 | 22.32927 | 923.4603 |
| C21 | 16.43286 | -4.155748 | -28.65575 | 22.46555 | 929.0966 |
| C21 | 16.45481 | -5.462670 | -29.96267 | 22.53091 | 931.7994 |
| C21 | 16.45315 | -5.363880 | -29.86388 | 22.52595 | 931.5943 |
| C21 | 16.44810 | -5.063625 | -29.56362 | 22.51090 | 930.9717 |
| C21 | 16.44502 | -4.880004 | -29.38000 | 22.50170 | 930.5916 |
| C21 | 16.44379 | -4.806987 | -29.30699 | 22.49805 | 930.4405 |
| C21 | 16.43948 | -4.550064 | -29.05006 | 22.48521 | 929.9096 |
| C21 | 16.45443 | -5.439936 | -29.93994 | 22.52977 | 931.7522 |
| C21 | 16.49080 | -7.598096 | -32.09810 | 22.63890 | 936.2654 |
| C21 | 16.45631 | -5.551668 | -30.05167 | 22.53538 | 931.9843 |
data_try <- expand.grid(
b = seq(170, 262, by = 10),
Ko = seq(2.1, 2.9, by = 0.2) * 10^-2,
eps_p = seq(16.5, 17.5, by = 1)
) %>% as_data_frame() %>%
mutate(
eps_f = case_when(eps_p ==16.5 ~ 24,
eps_p == 17.5 ~ 25),
compound = case_when (eps_p == 16.5 ~ "C17,19,21",
eps_p == 17.5 ~ "phytane"),
CO2_aq = b / (eps_f - eps_p),
pCO2 = CO2_aq / Ko
)
# base plot
base_plot <-
ggplot() +
aes(x = b, y = pCO2, linetype = factor(Ko)) + #removed color for eps_p because both eps_p lie on the same trend (becuase of eps_f assignment)
geom_line()
# plotting data 1
base_plot %+% data_try
cisotope <- read_excel(file.path("metadata", "Appendix_Table1_geochemistry.xlsx")) %>%
#rename columns
select(depth = `Abs. depth (m)` , d13c_org = `Average δ13Corg (‰ VPDB)` )
final_data <- full_join(cisotope, final_data, by = "depth")
timescale <- read_excel(file.path("metadata", "Appendix_Table2_SH1_agemodel.xlsx")) %>%
rename(depth = `Depth`) %>%
arrange(depth)
# extrapolation
extrapolated_timescale <-
data_frame(
depth = final_data$depth %>% unique()
) %>%
filter(!is.na(depth)) %>%
mutate(
in_range = depth <= max(timescale$depth) & depth >= min(timescale$depth),
Date.Ma = map_dbl(depth, function(d) {
fit_data <- mutate(timescale, closest = abs(d - depth)) %>% arrange(closest)
m <- lm(Date.Ma ~ depth, fit_data[1:10,])
p <- predict(m, data.frame(depth = d), se.fit = TRUE)
return(p$fit)
}) ,
in_range = depth <= max(timescale$depth) & depth >= min(timescale$depth),
ka = map_dbl(depth, function(d) {
fit_data <- mutate(timescale, closest = abs(d - depth)) %>% arrange(closest)
m <- lm(ka ~ depth, fit_data[1:10,])
p <- predict(m, data.frame(depth = d), se.fit = TRUE)
return(p$fit)
}))
ggplot() +
aes(depth, Date.Ma) +
geom_line(data = timescale) +
geom_point(data = extrapolated_timescale, mapping = aes(color = in_range))
alldata_w_time <- left_join(final_data, extrapolated_timescale, by = "depth")
alldata_w_time <- alldata_w_time %>%
mutate(age_range = case_when(.$compound %in% c("C17", "C19", "C21", "phytane") ~ (-0.00001) , #in MA (100 yrs)
.$compound %in% c("C27", "C29", "C31", "C33", "C35") ~ (-.01)), #in MA (10000 yrs)
age_dev = Date.Ma + age_range
)
alkane <- alldata_w_time %>%
filter(compound %in% c( "C29", "C31", "C33"), true_d13C_pred_in_range == TRUE, prep == "adduct")%>%
iso_plot_data(
x = Date.Ma,
y = true_d13C_pred,
#y_error = true_d13C_pred_se ,
shape = compound
#color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
scale_x_reverse(limits = c(94.5, 94.34)) +
#scale_x_continuous(limits = c(94.6, 94.3))+
#scale_y_continuous(limits = c(-33, -26)) +
coord_flip()
ggplotly(alkane)
bulkorg <- subset(alldata_w_time, ka > -100 & ka < 120) %>%
iso_plot_data(
x = ka,
y = d13c_org
) +
# additional features beyond the default plot
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
#scale_y_reverse() +
scale_x_continuous(limits = c(-100, 120)) +
#coord_cartesian(xlim=c(-20, 80))+
coord_flip()
grid.arrange(bulkorg, alkane, ncol = 2)
alkaneall <- alldata_w_time %>%
filter(compound %in% c("pristane", "phytane", "C19", "C21", "C31", "C33"), true_d13C_pred_in_range == TRUE)%>%
iso_plot_data(
x = Date.Ma,
y = true_d13C_pred,
y_error = true_d13C_pred_se ,
color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free", space = "free") +
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
scale_x_reverse() +
#scale_x_continuous()+
scale_y_continuous() +
coord_flip()
ggplotly(alkaneall)
bulkorgall <- subset(alldata_w_time) %>%
iso_plot_data(
x = ka,
y = d13c_org
) +
# additional features beyond the default plot
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
#scale_y_reverse() +
#scale_x_continuous(limits = c(-100, 120)) +
#coord_cartesian(xlim=c(-20, 80))+
coord_flip()
grid.arrange(bulkorgall, alkaneall, ncol = 2)
pools <- alldata_w_time %>%
filter(compound %in% c("C19", "C21", "C31", "C33"))%>%
ggplot()+
facet_grid(~compound, scales = "free", space = "free") +
geom_point(aes(x = Date.Ma , y = d13C_pool_avg)) +
#geom_errorbar(aes(x = ka , ymin = d13C_pool_min, ymax = d13C_pool_max)) +
#geom_errorbarh(aes(xmin = age_dev, xmax = ka , y = d13C_pool_avg)) +
#geom_line(aes(x = age_dev, y = d13C_pool_avg, color = "red")) +
#geom_line(aes(x = ka, y = d13C_pool_max)) +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "Inorganic pool d13C") +
scale_x_reverse() +
#scale_x_continuous()+
scale_y_continuous() +
#scale_x_continuous(limits = c(-100, 120)) +
#coord_cartesian(xlim=c(-20, 80))+
coord_flip()
ggplotly(pools)
pools_alttime <- alldata_w_time %>%
filter(Identifier.1 != "OG005_F1 alkanes_122m") %>% #removing duplicate for 94.45 Ma (non-adducted)
filter(compound %in% c( "C19", "C33"))%>%
ggplot()+
aes(x = Date.Ma , y = d13C_pool_avg, color = "red") +
facet_grid(~compound, scales = "free", space = "free") +
geom_line() +
#geom_errorbar(aes(x = ka , ymin = d13C_pool_min, ymax = d13C_pool_max)) +
#geom_errorbarh(aes(xmin = age_dev, xmax = ka , y = d13C_pool_avg)) +
#geom_line(aes(x = ka, y = d13C_pool_min)) +
geom_line(aes(x = age_dev, y = d13C_pool_avg, color = "blue")) +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "Inorganic pool d13C") +
scale_x_reverse(limits = c(94.6, 94.3)) +
#scale_x_continuous()+
scale_y_continuous() +
#scale_x_continuous(limits = c(94.6, 94.55)) +
#coord_cartesian(xlim=c(-20, 80))+
coord_flip()
ggplotly(pools_alttime)
alldata_w_time %>% filter(compound == "phytane") %>% select(Date.Ma, eps_p_lipid, d13C_pool_avg, true_d13C_pred, CO2_aq_conc, pCO2) %>% kable()
| Date.Ma | eps_p_lipid | d13C_pool_avg | true_d13C_pred | CO2_aq_conc | pCO2 |
|---|---|---|---|---|---|
| 94.45489 | 17.93042 | -6.505161 | -30.40516 | 24.04670 | 994.4872 |
| 94.30791 | 17.92757 | -6.349668 | -30.24967 | 24.03699 | 994.0855 |
| 94.34293 | 17.91180 | -5.490317 | -29.39032 | 23.98351 | 991.8741 |
| 94.22115 | 17.88922 | -4.257208 | -28.15721 | 23.90736 | 988.7246 |
| 94.37841 | 17.87707 | -3.592498 | -27.49250 | 23.86659 | 987.0385 |
| 94.17981 | 17.90731 | -5.245692 | -29.14569 | 23.96835 | 991.2470 |
| 94.04433 | 17.87767 | -3.625275 | -27.52528 | 23.86860 | 987.1214 |
| 94.56738 | 17.91410 | -5.615702 | -29.51570 | 23.99130 | 992.1959 |
| 94.52951 | 17.91765 | -5.809358 | -29.70936 | 24.00333 | 992.6935 |
| 93.63407 | 17.92682 | -6.309279 | -30.20928 | 24.03447 | 993.9813 |
| 93.91610 | 17.90231 | -4.972727 | -28.87273 | 23.95147 | 990.5487 |
| 93.66267 | 17.92183 | -6.037430 | -29.93743 | 24.01752 | 993.2804 |
| 94.27372 | 17.89758 | -4.714443 | -28.61444 | 23.93552 | 989.8892 |
| 93.85049 | 17.89374 | -4.504509 | -28.40451 | 23.92258 | 989.3540 |
| 94.10260 | 17.88144 | -3.831769 | -27.73177 | 23.88124 | 987.6445 |
| 94.01944 | 17.87633 | -3.551591 | -27.45159 | 23.86409 | 986.9350 |
| 93.78417 | 17.90732 | -5.246106 | -29.14611 | 23.96838 | 991.2481 |
| 93.73495 | 17.92271 | -6.085177 | -29.98518 | 24.02049 | 993.4034 |
| 93.68296 | 17.91013 | -5.399526 | -29.29953 | 23.97788 | 991.6412 |
| 93.75582 | 17.91882 | -5.872967 | -29.77297 | 24.00728 | 992.8571 |
| 93.61487 | 17.91623 | -5.731822 | -29.63182 | 23.99851 | 992.4942 |
| 93.57814 | 17.90877 | -5.324885 | -29.22488 | 23.97326 | 991.4499 |
| 93.81802 | 17.91766 | -5.810122 | -29.71012 | 24.00338 | 992.6954 |
| 94.13923 | 17.89399 | -4.517941 | -28.41794 | 23.92341 | 989.3882 |
#using my eps_p (lipd rather than caco3)
eps_p <- final_data %>%
filter(compound %in% c( "phytane" , "C19", "C21"), true_d13C_pred_in_range == TRUE)%>%
iso_plot_data(
x = depth,
y = eps_p_lipid
#y_error = true_d13C_pred_se ,
#shape = compound
#color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_point() +
#theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "eps_p") +
#scale_x_reverse() +
#scale_x_continuous(limits = c(-100, 120))+
#scale_y_continuous(limits = c(-33, -26)) +
coord_flip()
ggplotly(eps_p)
#using my eps_p (lipd rather than caco3)
pCO2 <- alldata_w_time %>%
filter(compound %in% c( "phytane" , "C19", "C21"), true_d13C_pred_in_range == TRUE)%>%
iso_plot_data(
x = Date.Ma,
y = pCO2
#y_error = true_d13C_pred_se ,
#shape = compound
#color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_point() +
#theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "pCO2") +
scale_x_reverse() +
#scale_x_continuous(limits = c(-100, 120))+
#scale_y_continuous(limits = c(-33, -26)) +
coord_flip()
ggplotly(pCO2)
pCO22 <- alldata_w_time %>%
filter(compound %in% c( "C21"), true_d13C_pred_in_range == TRUE)%>%
ggplot()+
aes(
x = Date.Ma,
y = pCO2
#y_error = true_d13C_pred_se ,
#shape = compound
#color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_smooth() +
geom_point() +
#theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "pCO2") +
scale_x_reverse(limits = c(94.5, 94.34)) +
#scale_x_continuous(limits = c(-100, 120))+
#scale_y_continuous(limits = c(-33, -26)) +
coord_flip()
alkane <- alldata_w_time %>%
filter(compound %in% c( "C31", "C33"), true_d13C_pred_in_range == TRUE, prep == "adduct")%>%
ggplot()+
aes(
x = Date.Ma,
y = d13C_pool_avg,
#y_error = true_d13C_pred_se ,
shape = compound
#color = prep
) +
# additional features beyond the default plot
facet_grid(~compound, scales = "free") +
geom_smooth() +
geom_point() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5)) +
labs(y = "d13C") +
scale_x_reverse(limits = c(94.5, 94.34)) +
#scale_x_continuous(limits = c(94.6, 94.3))+
#scale_y_continuous(limits = c(-33, -26)) +
coord_flip()
ggplotly(pCO2)
ggplotly(alkane)
#> `geom_smooth()` using method = 'loess' and formula 'y ~ x'
grid.arrange( alkane, pCO22, ncol = 2)
#> `geom_smooth()` using method = 'loess' and formula 'y ~ x'
#> `geom_smooth()` using method = 'loess' and formula 'y ~ x'
#> Warning: Removed 20 rows containing non-finite values (stat_smooth).
#> Warning: Removed 20 rows containing missing values (geom_point).